##################################################
## Author: Stephanie L. DeMora, Christian Lindke, Jenn Merolla, and Laura Stephenson
## Project: "Ready for a Woman President?" Paper Replication Script
## Date: September 25th 2020
##################################################

# Set up: Load in all packages
Packages <- c("ggplot2","hrbrthemes","tidyverse","haven","broom","Rmisc","dotwhisker","margins","jtools","mediation","here","rlist","dygraphs", "webshot", "htmlwidgets", "xts", "sensemakr")
lapply(Packages, library, character.only = TRUE)

# Plot Themes:
steph_theme <- list(
  hrbrthemes::theme_ipsum_rc(base_size = 12,grid= "",plot_title_size = 20,axis_text_size = 18),
  #ggsci::scale_color_material("grey"),
  #ggsci::scale_fill_material("grey"),
  theme(legend.position = "top", plot.caption = element_text(hjust = 0, face= "italic"),
        plot.title.position = "plot",
        plot.caption.position =  "plot"))

# Read in Data:
df <- read_dta("Data Files/CLEAN-CCES19_UCR.dta")
# Prep Data:
source("00. CCES_cleaning.R")


# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table A  --------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1$TREAT <- car::recode(df1$TREAT, "'Harder Win' = 'High Comfort 2'; 'Neighbor' = 'Low Comfort'; 'Personal' = 'High Comfort 1'")
options(scipen=999)
mod1 <- lm(UCR119b_1_1 ~ TREAT, data=df1)
mod2 <- lm(UCR119d_1_1 ~ TREAT, data=df1)
huxtable::huxreg(mod1,mod2)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table B  --------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1 %>%
  dplyr::select(UCR119b_1_1, TREAT) %>%
  dplyr::rename(DV = UCR119b_1_1) %>%
  mutate(Candidate = "Harris") -> dfSUB1

df1 %>%
  dplyr::select(UCR119d_1_1, TREAT) %>%
  dplyr::rename(DV = UCR119d_1_1) %>%
  mutate(Candidate = "Warren") -> dfSUB2

MASTER <- rbind(dfSUB1,dfSUB2)

Mod_EstDiff <- lm(DV ~ TREAT * Candidate, data=MASTER)
summary(Mod_EstDiff)


# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table C  --------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1$Democrat <- df1$party
df1$Woman <- df1$GENDER
mod1 <- lm(UCR119d_1_1 ~ TREAT * Woman * Democrat, data = df1)
mod2 <- lm(UCR119b_1_1 ~ TREAT * Woman * Democrat, data = df1)
export_summs(mod1,mod2, stars = c(`++` = 0.2, `*` = 0.1, `**` = 0.05),
               error_format = "({std.error})",
               error_pos = c("below"))


# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table D  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1$party <- as.numeric(as.character(labelled::remove_var_label(df1$party)))
mod1 <- lm(UCR119b_1_1 ~ TREAT * GENDER * party, data = df1)

margins_summary <- function(model, ..., level = 0.90, by_factor = TRUE) {
  mar <- margins(model, ...)
  summary(mar, level = 0.90, by_factor = by_factor)
}

tgc <- margins_summary(mod1, at = list(GENDER = range(df1$GENDER), party = range(df1$party)))
tgc$GENDER <- car::recode(tgc$GENDER, "0 = 'Men'; 1 = 'Women'")
tgc$party <- as.character(tgc$party)
tgc$model <- car::recode(tgc$party, "0 = 'Independent'; 1 = 'Democrat'")
tgc$term <- tgc$factor
tgc$conf.low = tgc$lower
tgc$conf.high = tgc$upper
tgc$conf.low_90 = tgc$lower
tgc$conf.high_90 = tgc$upper
tgc$estimate <- tgc$AME
margins_summary2 <- function(model, ..., level = 0.95, by_factor = TRUE) {
  mar <- margins(model, ...)
  summary(mar, level = 0.95, by_factor = by_factor)
}
tgc95 <- margins_summary2(mod1, at = list(GENDER = range(df1$GENDER), party = range(df1$party)))
tgc$conf.high_95 <- tgc95$upper
tgc$conf.low_95 <- tgc95$lower
tgc1 <- subset(tgc,tgc$GENDER == "Women")
tgc1$factor <- paste(tgc1$factor, " ", sep = "")
tgc1$term <- paste(tgc1$term, " ", sep = "")
tgc <- subset(tgc,tgc$GENDER == "Men")
master <- rbind(tgc,tgc1)
master %>%
  filter(str_detect(term, "TREAT")) -> master

master$term <- substring(master$term, 6)
master %>%
  arrange(GENDER) -> master

master$term <- car::recode(master$term, "'Harder Win' = 'High Comfort 2'; 'Neighbor' = 'Low Comfort'; 'Personal' = 'High Comfort 1'")
master$term <- car::recode(master$term, "'Harder Win ' = 'High Comfort 2 '; 'Neighbor ' = 'Low Comfort '; 'Personal ' = 'High Comfort 1 '")
master$sort <- ""
master$sort[master$term == "Low Comfort"] <- 1
master$sort[master$term == "Low Comfort "] <- 4
master$sort[master$term == "High Comfort 1"] <- 2
master$sort[master$term == "High Comfort 1 "] <- 5
master$sort[master$term == "High Comfort 2"] <- 3
master$sort[master$term == "High Comfort 2 "] <- 6

master %>%
  arrange(sort) -> master
master$term <- factor(master$term,levels = c("High Comfort 2 ","High Comfort 1 ","Low Comfort ","High Comfort 2","High Comfort 1","Low Comfort"))
master %>%
  dplyr::select(term,GENDER,party,AME,SE,z,p) -> For_Table_DEMIND
For_Table_DEMIND$party <- car::recode(For_Table_DEMIND$party, "0 = 'Independent'; 1 = 'Democrat'")

# For Republicans now:
df %>%
  filter(PARTISANSHIP != "Democrat") %>%
  filter(!is.na(GENDER))-> df1_new
df1_new$TREAT <- car::recode(df1_new$TREAT, "'Harder Win' = 'High Comfort 2'; 'Neighbor' = 'Low Comfort'; 'Personal' = 'High Comfort 1'")
df1_new$PARTISANSHIP2 <- as.numeric(as.character(car::recode(df1_new$PARTISANSHIP, "'Independent' = 0; 'Republican' = 1")))
linked30 <- lm(UCR119b_1_1 ~ TREAT*GENDER*PARTISANSHIP2, data=df1_new, na.action = na.omit)
rep_ind <- summary(margins(linked30, at = list(GENDER = range(df1_new$GENDER), PARTISANSHIP2 = range(df1_new$PARTISANSHIP2))))
rep_ind %>%
  filter(str_detect(factor, "TREAT")) -> rep_ind
rep_ind$factor <- substring(rep_ind$factor, 6)
rep_ind %>%
  arrange(GENDER) %>%
  filter(PARTISANSHIP2 == 1) %>%
  mutate(PARTISANSHIP2 = "Republican") %>%
  dplyr::rename(term = factor,
                party = PARTISANSHIP2) -> rep_ind
rep_ind$GENDER <- car::recode(rep_ind$GENDER, "0 = 'Men'; 1 = 'Women'")
rep_ind %>%
  dplyr::select(term,GENDER,party,AME,SE,z,p) -> For_Table_REP

Table_D <- rbind(For_Table_REP,For_Table_DEMIND)
Table_D

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table E  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1$party <- as.numeric(as.character(labelled::remove_var_label(df1$party)))
mod1 <- lm(UCR119d_1_1 ~ TREAT * GENDER * party, data = df1)

margins_summary <- function(model, ..., level = 0.90, by_factor = TRUE) {
  mar <- margins(model, ...)
  summary(mar, level = 0.90, by_factor = by_factor)
}

tgc <- margins_summary(mod1, at = list(GENDER = range(df1$GENDER), party = range(df1$party)))
tgc$GENDER <- car::recode(tgc$GENDER, "0 = 'Men'; 1 = 'Women'")
tgc$party <- as.character(tgc$party)
tgc$model <- car::recode(tgc$party, "0 = 'Independent'; 1 = 'Democrat'")
tgc$term <- tgc$factor
tgc$conf.low = tgc$lower
tgc$conf.high = tgc$upper
tgc$conf.low_90 = tgc$lower
tgc$conf.high_90 = tgc$upper
tgc$estimate <- tgc$AME
margins_summary2 <- function(model, ..., level = 0.95, by_factor = TRUE) {
  mar <- margins(model, ...)
  summary(mar, level = 0.95, by_factor = by_factor)
}
tgc95 <- margins_summary2(mod1, at = list(GENDER = range(df1$GENDER), party = range(df1$party)))
tgc$conf.high_95 <- tgc95$upper
tgc$conf.low_95 <- tgc95$lower
tgc1 <- subset(tgc,tgc$GENDER == "Women")
tgc1$factor <- paste(tgc1$factor, " ", sep = "")
tgc1$term <- paste(tgc1$term, " ", sep = "")
tgc <- subset(tgc,tgc$GENDER == "Men")
master <- rbind(tgc,tgc1)
master %>%
  filter(str_detect(term, "TREAT")) -> master

master$term <- substring(master$term, 6)
master %>%
  arrange(GENDER) -> master

master$term <- car::recode(master$term, "'Harder Win' = 'High Comfort 2'; 'Neighbor' = 'Low Comfort'; 'Personal' = 'High Comfort 1'")
master$term <- car::recode(master$term, "'Harder Win ' = 'High Comfort 2 '; 'Neighbor ' = 'Low Comfort '; 'Personal ' = 'High Comfort 1 '")
master$sort <- ""
master$sort[master$term == "Low Comfort"] <- 1
master$sort[master$term == "Low Comfort "] <- 4
master$sort[master$term == "High Comfort 1"] <- 2
master$sort[master$term == "High Comfort 1 "] <- 5
master$sort[master$term == "High Comfort 2"] <- 3
master$sort[master$term == "High Comfort 2 "] <- 6

master %>%
  arrange(sort) -> master
master$term <- factor(master$term,levels = c("High Comfort 2 ","High Comfort 1 ","Low Comfort ","High Comfort 2","High Comfort 1","Low Comfort"))
master %>%
  dplyr::select(term,GENDER,party,AME,SE,z,p) -> For_Table_DEMIND
For_Table_DEMIND$party <- car::recode(For_Table_DEMIND$party, "0 = 'Independent'; 1 = 'Democrat'")

# For Republicans now:
df %>%
  filter(PARTISANSHIP != "Democrat") %>%
  filter(!is.na(GENDER))-> df1_new
df1_new$TREAT <- car::recode(df1_new$TREAT, "'Harder Win' = 'High Comfort 2'; 'Neighbor' = 'Low Comfort'; 'Personal' = 'High Comfort 1'")
df1_new$PARTISANSHIP2 <- as.numeric(as.character(car::recode(df1_new$PARTISANSHIP, "'Independent' = 0; 'Republican' = 1")))
linked30 <- lm(UCR119d_1_1 ~ TREAT*GENDER*PARTISANSHIP2, data=df1_new, na.action = na.omit)
rep_ind <- summary(margins(linked30, at = list(GENDER = range(df1_new$GENDER), PARTISANSHIP2 = range(df1_new$PARTISANSHIP2))))
rep_ind %>%
  filter(str_detect(factor, "TREAT")) -> rep_ind
rep_ind$factor <- substring(rep_ind$factor, 6)
rep_ind %>%
  arrange(GENDER) %>%
  filter(PARTISANSHIP2 == 1) %>%
  mutate(PARTISANSHIP2 = "Republican") %>%
  dplyr::rename(term = factor,
                party = PARTISANSHIP2) -> rep_ind
rep_ind$GENDER <- car::recode(rep_ind$GENDER, "0 = 'Men'; 1 = 'Women'")
rep_ind %>%
  dplyr::select(term,GENDER,party,AME,SE,z,p) -> For_Table_REP

Table_E <- rbind(For_Table_REP,For_Table_DEMIND)
Table_E

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table F  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1 %>%
  filter(UCR123_1 > 1) %>%
  filter(UCR123_1 != 8) -> df1_likely

df1_likely$party <- as.numeric(as.character(labelled::remove_var_label(df1_likely$party)))

mod1 <- lm(UCR119b_1_1 ~ TREAT * GENDER * party, data = df1_likely)
summary(mod1)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table G  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
mod1 <- lm(UCR119d_1_1 ~ TREAT * GENDER * party, data = df1_likely)
summary(mod1)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table H  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
mod1 <- lm(UCR119b_1_1 ~ TREAT * GENDER * party, data = df1_likely)
tgc <- margins_summary(mod1, at = list(GENDER = range(df1_likely$GENDER), party = range(df1_likely$party)))
tgc[-c(1:8),]

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table I  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
mod1 <- lm(UCR119d_1_1 ~ TREAT * GENDER * party, data = df1_likely)
tgc <- margins_summary(mod1, at = list(GENDER = range(df1_likely$GENDER), party = range(df1_likely$party)))
tgc[-c(1:8),]


# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table J  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1 %>%
  filter(race != 1) -> POC
# Harris
mod1 <- lm(UCR119b_1_1 ~ TREAT, data = POC, subset = party != 99)
# Warren:
mod2 <- lm(UCR119d_1_1 ~ TREAT, data = POC, subset = party != 99)

export_summs(mod1,mod2, model.names = c("Harris","Warren"), ci_level = 0.9)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table K  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
# Harris - Low Comfort - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_MEN_DEMS.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - Low Comfort - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_MEN_IND.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - Low Comfort - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_WOMEN_DEMS.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - Low Comfort - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_WOMEN_IND.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_MEN_DEMS_HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_MEN_IND-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_WOMEN_DEMS-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_WOMEN_IND-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_MEN_DEMS_HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_MEN_IND-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_WOMEN_DEMS-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_HarrisAlone_WOMEN_IND-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_MEN_DEM.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_MEN_IND.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_WOMEN_DEM.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_WOMEN_IND.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_MEN_DEM-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_MEN_IND-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_WOMEN_DEM-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_WOMEN_IND-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_MEN_DEM-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_MEN_IND-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_WOMEN_DEM-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1b_Results_WarrenAlone_WOMEN_IND-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table L  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
# Harris - Low Comfort - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Dems_MEN.RData")
#plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - Low Comfort - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Ind_MEN.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - Low Comfort - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Dems_WOMEN.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - Low Comfort - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Ind_WOMEN.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Dems_MEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Ind_MEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Dems_WOMEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 1 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Ind_WOMEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Dems_MEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Ind_MEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Dems_WOMEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Harris - High Comfort 2 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed1_Results_Ind_WOMEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Dems_MEN.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Ind_MEN.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Dems_WOMEN.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - Low Comfort - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Ind_WOMEN.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Dems_MEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Independent Men
results_modmed1 <- list.load("FigureModMed2_Results_Ind_MEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Dems_WOMEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 1 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Ind_WOMEN-HC1.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Democratic Men
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Dems_MEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Independent Men
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Ind_MEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Democratic Women
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Dems_WOMEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# Warren - High Comfort 2 - Independent Women
results_modmed1 <- list.load("Data Files/FigureModMed2_Results_Ind_WOMEN-HC2.RData")
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9)
summary(results_modmed1)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table M  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
mod1 <- lm(UCR119a_1_1 ~ TREAT, data = df1)
mod2 <- lm(UCR119c_1_1 ~ TREAT, data = df1)
export_summs(mod1,mod2, model.names = c("Biden","Sanders"), ci_level = 0.9)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table N  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
mod1 <- lm(UCR119a_1_1 ~ TREAT * GENDER * party, data = df1)
mod2 <- lm(UCR119c_1_1 ~ TREAT * GENDER * party, data = df1)
export_summs(mod1,mod2, model.names = c("Biden","Sanders"), ci_level = 0.9)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table O  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
mod1 <- lm(UCR119a_1_1 ~ TREAT * GENDER * party, data = df1)
tgc <- margins_summary(mod1, at = list(GENDER = range(df1$GENDER), party = range(df1$party)))
tgc[-c(1:8),]

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table P  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
mod1 <- lm(UCR119c_1_1 ~ TREAT * GENDER * party, data = df1)
tgc <- margins_summary(mod1, at = list(GENDER = range(df1$GENDER), party = range(df1$party)))
tgc[-c(1:8),]

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Table Q  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
df1$TREAT2 <- relevel(df1$TREAT,ref = "High Comfort 1")
# Harris
data2 <- tidy(lm(UCR119b_1_1 ~ TREAT2, data=df1), conf.level = 0.90, conf.int = TRUE)
data2$p.value <- as.character(round(data2$p.value,digits = 3))
data2[3,]

# Warren:
data4 <- tidy(lm(UCR119d_1_1 ~ TREAT2, data=df1), conf.level = 0.90, conf.int = TRUE)
data4$p.value <- as.character(round(data4$p.value,digits = 3))
data4[3,]

# BIDEN
data4 <- tidy(lm(UCR119a_1_1 ~ TREAT2, data=df1), conf.level = 0.90, conf.int = TRUE)
data4$p.value <- as.character(round(data4$p.value,digits = 3))
data4[3,]

# SANDERS
data4 <- tidy(lm(UCR119c_1_1 ~ TREAT2, data=df1), conf.level = 0.90, conf.int = TRUE)
data4$p.value <- as.character(round(data4$p.value,digits = 3))
data4[3,]

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure A  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
data <- data.frame(
  time=seq(from=Sys.Date()-582, to=Sys.Date()-554, by=1 ), 
  Warren=c(21.3,21.3,21.3,21.3,21.3,20.7,20.7,20.9,20.9,20.9,20.5,20.5,20.2,19.8,19.7,19.7,17.8,17.2,17.6,17.1,16.6,16.4,16.4,16.4,16.0,16.0,15.7,15.7,15.5),
  Harris=c(4.5,4.5,4.5,4.5,4.5,4.5,4.5,4.5,4.5,4.5,4.5,4.5,4.4,4.5,4.5,4.2,4.1,4.1,4.1,4.0,3.9,3.9,3.9,3.9,4.0,4.0,3.9,3.9,3.9),
  Biden=c(26.9,26.9,26.9,26.9,26.9,26.9,26.9,26.7,26.6,26.6,26.3,26.3,26.4,26.3,26.7,26.7,26.8,26.7,26.7,26.6,26.3,26.2,26.2,26.2,26.0,26.0,26.0,26.0,26.0),
  Sanders=c(15.8,15.8,15.8,15.8,15.8,16.0,16.0,16.0,15.9,15.9,16.2,16.2,16.3,16.4,16.1,16.8,17.2,17.2,17.2,17.2,17.0,16.9,16.9,16.9,16.9,16.9,16.8,16.8,15.5))
don=xts( x=data[,-1], order.by=data$time)
p <- dygraph(don,main = "Primary Candidate National Polling <br> <small>FiveThirtyEight Data</small>",
             ylab = "Percent Averages",
             xlab = "Study 1 Field Dates") %>%
  dyLegend(width = 135)
p

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure B1  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
results_modmed1 <- list.load("Data Files/FigureMod_Results_HarrisAlone-HC1.RData")
summary(results_modmed1)
plot(results_modmed1,sub = "90% Confidence Intervals")

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure B2  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
results_modmed1 <- list.load("Data Files/FigureMod_Results_HarrisAlone-HC2.RData")
summary(results_modmed1)
plot(results_modmed1,sub = "90% Confidence Intervals")

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure B3  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
results_modmed2 <- list.load("Data Files/FigureMod_Results_WarrenAlone-H1.RData")
summary(results_modmed2)
plot(results_modmed2,sub = "90% Confidence Intervals")

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure B4  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
results_modmed2 <- list.load("Data Files/FigureMod_Results_WarrenAlone-H2.RData")
summary(results_modmed2)
plot(results_modmed2,sub = "90% Confidence Intervals")

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure C1  -------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
results_modmed1 <- list.load("Data Files/Figure4_Results.RData")
results_modmed2 <- list.load("Data Files/Figure5_Results.RData")

par(mfrow = c(1,2))
plot(results_modmed1,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9, adj = .5, main = "Panel A: Harris")
plot(results_modmed2,sub = "90% Confidence Intervals", cex.sub = 1, cex.axis = 1.27, adj = .9, adj = .5, main = "Panel B: Warren")
dev.off()

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure C2  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
set.seed(12082020)
df2 <- df[!is.na(df$Gender_Vote),]
df2$TREAT <- car::recode(df2$TREAT, "'Harder Win' = 'High Comfort 2'; 'Neighbor' = 'Low Comfort'; 'Personal' = 'High Comfort 1'")
model.0 <- lm(Gender_Vote ~ TREAT, df2)
model.M <- lm(UCR119b_1_1 ~ TREAT, df2)
model.Y <- lm(Gender_Vote ~ TREAT + UCR119b_1_1, df2)
results <- mediate(model.M, model.Y, treat='TREAT', mediator='UCR119b_1_1',
                   boot=TRUE, sims=1000, conf.level = .90, treat.value = "'High Comfort 1'", control.value = "Control")

plot(results)
summary(results)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure C3  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
set.seed(12082020)
df2 <- df[!is.na(df$Gender_Vote),]
df2$TREAT <- car::recode(df2$TREAT, "'Harder Win' = 'High Comfort 2'; 'Neighbor' = 'Low Comfort'; 'Personal' = 'High Comfort 1'")

model.0 <- lm(Gender_Vote ~ TREAT, df2)
model.M <- lm(UCR119b_1_1 ~ TREAT, df2)
model.Y <- lm(Gender_Vote ~ TREAT + UCR119b_1_1, df2)
results <- mediate(model.M, model.Y, treat='TREAT', mediator='UCR119b_1_1',
                   boot=TRUE, sims=1000, conf.level = .90, treat.value = "High Comfort 2", control.value = "Control")

plot(results)
summary(results)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure C4  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
model.0 <- lm(Gender_Vote ~ TREAT, df2)
model.M <- lm(UCR119d_1_1 ~ TREAT, df2)
model.Y <- lm(Gender_Vote ~ TREAT + UCR119d_1_1, df2)
results <- mediate(model.M, model.Y, treat='TREAT', mediator='UCR119d_1_1',
                   boot=TRUE, sims=1000, conf.level = .90, treat.value = "High Comfort 1", control.value = "Control")

plot(results)
summary(results)

# ------------------------------------------------------------------------------------------------------------------------
# ----------------------------------------------  Figure C5  ---------------------------------------------------------------
# ------------------------------------------------------------------------------------------------------------------------
model.0 <- lm(Gender_Vote ~ TREAT, df2)
model.M <- lm(UCR119d_1_1 ~ TREAT, df2)
model.Y <- lm(Gender_Vote ~ TREAT + UCR119d_1_1, df2)
results <- mediate(model.M, model.Y, treat='TREAT', mediator='UCR119d_1_1',
                   boot=TRUE, sims=1000, conf.level = .90, treat.value = "High Comfort 2", control.value = "Control")

plot(results)
summary(results)


